Neel and Valence Bond Crystal Order on a Distorted Kagome Lattice: Implications 

For Zn-Paratacamite 
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Zn-Paratacamite is a rare spin 1/2 antiferromagnetic insulator with an ideal kagome lattice struc- 
ture in part of its phase diagram. As a function of Zn doping, this material undergoes a structural 
distortion which relieves the frustration and introduces magnetic order in the ground state, though 
the precise nature of the order is not clear at this point. In this paper, we present strong evidence 
for Neel ordering in the strongly distorted phase of Zn-Paratacamite through the apphcation of 
quantum Monte-Carlo techniques. These numerical results support a recent Schwinger-boson mean 
field theory of Zn-Paratacamite. For weak distortion, close to the ideal kagome limit, our results 
indicate a regime with no Neel order but with a broken glide-plane symmetry. For this model the 
glide-plane symmetry is broken by any valence bond crystal. Hence, our results lend support to 
recent proposals^'? of a valence bond crystal ground state for the undistorted lattice. The phase 
transition between the two phases could be in the deconfined universality class if it is not a first 
order transition. 

PACS numbers: 75.10.Jm, 73.43.Nq, 75.40.Mg 



I. INTRODUCTION 



Recently, a number of spin 1/2 frustrated magnetic 
insulators have been discovered without any sign of mag- 
netic order or structural distortions down to the low- 
est temperatures studied^i^i^i^. Among these materials, 
the Zn doped Paratacamite family stands out for having 
a (nearly) controllable degree of distortion allowing the 
amount of geometric frustration to be tuned directly by 
an experimentalist. As such, they are a promising place 
to look for new phases of matter while at the same time 
probe how these new phases may be related to more well 
understood phases. 

The control of the distortion is largely through the 
chemical pressure induced by the substitution of Zn 
atoms for Cu atoms on the (grey) triangular lattice planes 
that live in between kagome-planes, as shown in Fig. [1] 
While Zn and Cu atoms are similar in size, Zn atoms 
fit into these sites without disrupting their environment, 
unlike Cu atoms which distort the kagome plans given a 
high enough density. In particular, for less than 0.3 fill- 
ing of Zn atoms (greater than 0.7 filling of Cu atoms) the 
lattice distorts in a remarkable bi-partite structure and 
magnetic order is found in the ground state^i^i^. The 
spins are thus relatively unfrustrated at these low doping 
concentrations. For Zn doping larger than this threshold, 
the lattice has the undistorted ideal kagome form and for 
X > 0.4 no magnetic ordering has been reported down to 
50 mK despite an estimated spin exchange J ~ 200K—. 

A natural theoretical model of this material is the spin 
1/2 Heisenberg model with two exchange parameters on 




FIG. 1: (Color online) The layered Zn Paratacamite magnetic 
lattice structure. Cu atoms live on the (sometimes distorted) 
kagome layer (bronze atoms) while Zn or Cu atoms occupy 
sites on a triangular lattice above the kagome plane (grey 
atoms). Solid (purple) bonds represent nearest neighbors {{ij) 
in Eq. ((T))) resulting in a 'brick-wall' lattice, while transparent 
(bronze) bonds represent the next nearest neighbors {{{ij}) in 
Eq. (HI). 



a distorted kagome lattice^i^ii^ (see Fig. [T|): 



(1) 



Here Ad tunes the distortion on the next nearest neighbor 
bonds (nnn) with Xd — I the undistorted ideal kagome 
limit (where nnn bonds are equivalent to nn bonds) . This 
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is an idealized model for several reasons: it assumes the 
coupling between planes and further neighbors is weak 
(which seems reasonable^), it neglects Dzyaloshinsky- 
Moria interactions possibly important for the low tem- 
perature susceptibility^'^: , and it replaces the effect of 
doping in the triangle lattice planes with the uniform dis- 
tortion parameter Ad. While it is possible that any of 
these approximations may be important for some prop- 
erties of Zn Paratacamite, here we will focus on those 
properties which clearly belong to the phenomenology of 
this simplified model. 

In this paper, we study the ground state properties of 
the Hamiltonian in Eq. ([T]) as a function of Ad , extrapo- 
lating between a bi-partite 'brick- wall' lattice at Ad = 
and the isotropic kagome lattice at Ad = 1- At Ad = 0, 
we show using valence bond quantum Monte-Carlo ^^'^'^ , 
that the ground state is magnetically ordered with the 
expected Neel pattern for this bi-partite lattice and with 
a magnetization of = 0.240(1) that is 22% smaller 
than the square lattice valu&i^. For < Ad < 1 we 
study this model using exact diagonalization on finite 
size clusters of size 12, 24 and 36 sites. By introducing 
symmetry breaking fields, we study the susceptibility of 
the ground state towards dimerization. Remarkably, we 
find that for Ad ^ 0.8, a phase transition occurs towards 
a rotationally invariant state which prefers to have a bro- 
ken glide-plane symmetry, consistent with the presence of 
a VBC order including the pin-wheel VBC pattern pro- 
posed by Ref. This symmetry breaking survives up 
to the Ad = 1 ideal kagome limit. While it is difficult 
to draw definitive conclusions on such small systems, a 
broken glide-plane symmetry support Refs. mm pro- 
posal that the spin 1/2 kagome antiferromagnet has a va- 
lence bond crystal (VBC) ground state. At the same time 
a broken glide-plane symmetry is not consistent with a 
spin-liquid phase, frequently supported by other exact di- 
agonalization studies^^. In addition, while we can't rule 
out a first order transition from a VBC phase to the Neel 
phase in our model, it is also possible this quantum phase 
transition is in the deconfined universality class^^. 



II. RESULTS AT Ad = 

We first discuss our results obtained at Ad = where 
we have been able to study large system. As can be 
seen from Figs. [U where the transparent (bronze) bonds 
are proportional to Ad, the lattice formed by the re- 
maining solid (purple) bonds is a bi-partite 'brick-wall' 
lattice with a coordination number of 3 on two thirds 
of the sites and of 2 on the remaining sites. Due to 
the bi-partite nature of the Ad = lattice there is no 
frustration. A classical AF Neel state can be unam- 
biguously assigned to the lattice. It is then possible to 
perform very efficient quantum Monte Carlo simulations 
using the recently proposedi2,il2, valence bond quantum 
Monte Carlo (VBQMC). For Ad = there is no sign 
problem and extremely precise results can be obtained. 



0.040 1 

0.035 r ^ 
0.030 



0.025 r 




1/N 

FIG. 2: (Color online.) The singlet-triplet gap, A, versus 
inverse system size 1/A''. The different curves correspond to 
different expansion orders n, with the ratio n/Ni, kept fixed. 
Here, Nb is the number of bond operators in the Hamiltonian, 
Ni = (4/3)iV. Results are shown for n/Nt = 5 (■), n/Nt = 7 
(o), n/Nb = 10 (•). Error bars are shown but are typically 
smaller than the symbols. 



VBQMC is a projection method where the T = ground- 
state is projected out through the repeated application 
of the hamiltonian, H, on a trial state, |4't). In essence, 
I'^g) = (-^f)"|*T)- In the limit where n ^ oo this 
becomes exact. In a practical implementation n is kept 
fixed at a high number and the different terms in |5'g) are 
sampled using Monte Carlo methods. For convergence, 
the relevant lattice size independent expansion order is 
n/Nh where TV;, is the number of terms in the Hamilto- 
nian. Nt is equal to (4/3)7V for the brick-wall lattice with 
N the number of sites in the lattice. Typically we use 
n/Nb = 3—10 and an extrapolation to n/Nb = oo can 
then be performed. 

We have performed systematic VBQMC studies of 
Ad = brick- wall lattices with number of sites N — 
12 X rn^ for to = 1, ... 10 using periodic boundary con- 
ditions. Typically, 10^ — 10^ MCS were performed for a 
range of values of n/Nb = 3, 5, 7, 10. All errrorbars were 
calculated using standard binning techniques. 

A very natural question to ask is if the Ad = brick- 
wall lattice has a non-zero singlet to triplet gap, A. A 
particularly appealing feature of VBQMC is that it al- 
lows for a direct estimator-i^ii^ of this gap independent 
of the estimators for the ground-state singlet and ex- 
cited triplet energies. Due to a cancellation of errors 
it is then possible to calculate this gap with a precision 
significantly exceeding that which could have been ob- 
tained by separately calculating the ground and excited 
state energies. Our results for A at Ad = are shown in 
Fig. [21 Data arc shown for 3 different values of n/Nb — 5 
(■), n/Nb = 7 (o), n/Nb = 10 (•) versus inverse sys- 
tem size 1/iV. At TV = 1200, n/Nb = 10 the gap is 
A ~ 0.0059(1) J and minimal dependence on the expan- 
sion order n/Ns is seen. From the results shown in Fig. [2] 
we conclude that the gap vanishes in the thermodynamic 
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FIG. 3: (Color online.) The structure factor, 5(qc), ver- 
sus inverse linear system size l/\/N. Results are shown for 
n/Nb = 10. Error bars are shown but are typically smaller 
than the symbols. Results for n/Ni, = 7 are indistinguishable 
from the n/Nb ~ 10 results shown and have been left out for 
clarity. 



limit. 

For the two-dimensional square lattice anti- 
ferromagnetic Heisenberg model it is well knowni^ 
that the anti-ferromagnetic order exists at T = with 
= 0.30743(1)^'*. The square lattice has a coordina- 
tion number of 4 where as the brick-wall lattice has a 
mixed coordination of 2 and 3. We therefore expect 
to be smaller or possibly zero for the brick-wall lattice. 
As usual, we define: 



S{qc 



\ \ x,y / 



(2) 



where Qc is the wave-vector of the staggered magnetiza- 
tion and S^(x, y) is given by: 



(3) 



with ex,y equal to -f 1 or — 1 depending on what sublattice 
the point (x, y) belongs to. Hence we havei^: 



(St) = lim V3S{^) 



(4) 



Our results for 5(qc) for n/Ni, — 10 are shown in Fig. [31 
It is expected^* that the leading finite size corrections are 
of the form and a fit to this form yields 5'(qc) = 

0.0192(5) and consequently: 



ml' = 0.240(3). 



(5) 



As expected, this value is reduced with respect to the 
square lattice result, but is clearly non-zero, indicating a 
well established anti- ferromagnetic order at A^i = 0. 




FIG. 4: (Color online.) (a) The glide-plane symmetry break- 
ing field, (b) C2 symmetry breaking field. Solid dimers denote 
bonds where the coupling strength is modified to J'. 



III. RESULTS AT Ad / 

We now turn to a discussion of our results for < 
Aij < 1. In this case it is no longer possible to perform 
VBQMC calculations due to a sign problem that appears 
rather severe as soon as Xjj ^ and reliable numerical 
results are therefore much harder to obtain. In light of 
the strong sign problem we have performed exact diag- 
onalization studies for < A < 1 on finite size systems 
employing periodic boundary conditions. Our goal is to 
study generalized bond susceptibilities with respect to 
symmetry breaking fields. We focus on C2 and glide- 
plane (GP) symmetry breaking fields shown in Fig. [4] 
where the dimers indicate bonds where the strength is 
modified J' ~ J ± S. The C2 symmetry corresponds to 
a rotation by tt and clearly the field shown in Fig. [H^a) 
breaks this symmetry. The GP symmetry^ is somewhat 
more exotic and corresponds to a translation along the 
rails where the dimers are sitting followed by a reflection 
around one of these rails. We note that the GP field does 
not break the C2 symmetry and likewise the C2 field pre- 
serves the GP symmetry. The pin-wheel VBC discussed 
in Ref. 9 would break the GP symmetry but not the 
C2 symmetry where as the columnar VBC^ would break 
both. If we by b,bc2,GP denote the ground-state ex- 
pectation value {Si ■ Sj) for the bond b and its partner 
under the symmetry operation bcj.GPj we can define the 
generalized bond susceptibility as follows: 

XC2,GP = 

|Abc„Gp(J' = J + 5) - Abc„Gp(J' ^J~S)\ 



lim 

5^0 



with 
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Abc,,GP =b(J')-bc,,Gp("f)■ 



(6) 



(7) 



Clearly, if Ab goes to zero linearly with 6 the generalized 
bond susceptibility is a constant independent of system 
size and the associated symmetry is not spontaneously 



FIG. 5: (Color online.) The 12, 24 and 36 site lattices. 



broken. On the other hand, a bond susceptibility diverg- 
ing with system size signals that the associated symmetry 
is spontaneously broken in the thermodynamic limit. 

When performing exact diagonalization studies of 
small systems the choice of the finite cluster is crucial 
since the smaller clusters will reduce the point group 
symmetry of the infinite lattice. For the isotropic kagome 
lattice, Xd — I, the plane group is p6mm. This symme- 
try group implies that for the isotropic kagome lattice 
all bonds are equivalent. Our choice of finite clusters are 
shown in Fig. E] for = 12, 24, 36. Only the = 36 
cluster has the full symmetry point group symmetry of 
the kagome lattice. However, for all clusters do we find 
that all bonds are equivalent at Xjj = 1. For < A^i < 1 
only two different type of bonds occur for these clusters. 
These are the only clusters we have found with these 
properties. For the bond susceptibilities to yield mean- 
ingful information about the thermodynamic limit this 
is very important since we want to make sure that the 
presence of a reduced point group symmetry doesn't ex- 
plicitly break the C2 or GP symmetry. This is not the 
case for the clusters shown in Fig. O 

Using the clusters from Fig. [5] we can now study xgp 
and XC2 a function of Ad for the different clusters. Our 
results are shown in Figs I6I7I respectively. We have used 
S < 0.001 small enough that x, Eq. is almost com- 
pletely independent of S. We begin by discussing xgp 
shown in Fig. [6] For Ad less than roughly ^ 0.8 do we 
find that xgp is almost independent of N. In the inset 
is shown 1/xgp as a function of 1/A^ for A — 0.3 indi- 
cating a finite value in the thermodynamic limit. This 
is consistent with the GP symmetry not being broken. 
However, for Ad greater than ~ 0.8 pronounced size de- 
pendence occurs. At Ad = 1, 1/xgp as a function of 
1/A^ is shown in the inset. In this case it seems reason- 
able to conclude that xgp diverges with N and hence 
that the GP symmetry is spontaneously broken in the 
thermodynamic limit. A natural interpretation of this 
result is that a quantum phase transition occurs at 0.8 



FIG. 6: (Color online.) The glide-plane susceptibility, xgp as 
a function of Ao for the different system sizes, A'^ — 12, 24, 36. 
The circles represent results for A'' = 36, A_d = 0.3. The inset 
shows Xgp versus 1/N at Ad = 0.3, 1.0. 




FIG. 7: (Color online.) The C2 susceptibility, XC2 as a func- 
tion of Ad for the different system sizes. A' = 12, 24, 36. The 
circles represent results for A' — 36, Ad = 0.3. The inset 
shows XC2 versus 1/N at Ad = 0.3, 1.0. 



between a state with anti-ferromagnetic order, that does 
not break the GP symmetry, to a new phase where the 
GP symmetry is broken. 

Finally, in Fig.[7]we show our results for XC2 ■ Again we 
see that for Ad smaller than roughly ~ 0.8 there is very 
little A^ dependence. In the inset in Fig.[7]is shown I/XC2 
as a function of 1/A^ at Ad = 0.3. Clearly the results ex- 
trapolate to a finite value in the thermodynamic limit 
consistent with the absence of C'2 symmetry breaking as 
would be the case for an anti-ferromagnetic phase. As be- 
fore, we find that for Ad greater than roughly ~ 0.8 pro- 
nounced finite size effects occur consistent with a quan- 
tum phase transition. However, in this case, as can be 
seen in the inset in Fig. [7] at Ad = 1, the susceptibility 
doesn't diverge but rather tends to a finite, possibly very 



5 



m 

(a)MZ pattern (b)Pin-wheel pattern 

FIG. 8: (a) The Marston and Zeng(MZ) dimer patternia on 
the kagome lattice and (b) the pin-wheel state on the dis- 
torted kagome lattice^.. The MZ pattern arises by maximiz- 
ing the number of dimers around each hexagon. Notice how 
the dimers on the two highlighted hexagons can be rotated by 
180 about the center of the hexagon without needing to alter 
the rest of the dimer pattern. This benzene-like resonance 
suggests that rotational symmetry may be only very weakly 
broken if this is the ground state. The pin-wheel pattern, on 
the other hand, maximizes dimers around each rhombus and 
is manifestly C2 rotationally symmetric about each hexagon. 

small value in the thermodynamic limit. Therefore, in 
the new phase occurring for Xd greater than 0.8 the C2 
symmetry is not broken. 

IV. DISCUSSION 

Given these results, we may draw several conclusions. 
At \d ~ 0, we have demonstrated that the spin gap 
vanishes and that the ground state has a finite staggered 
magnetization that is 22% smaller than the square lattice 
value. We find this surprising given that each site has 
either two or three neighbors (with 2.67 neighbors on 
average) and that this network is not far from the one 
dimensional chain model which is disordered. It appears 
that this form of dimensional reduction does not easily 
suppress magnetic order. 

While we could only study very small system sizes for 
Ad > 0, finite size effects seem to be small all the way out 
to Ajj w 0.8. As a result, the anti-ferromagnetic order is 
quite robust and appears to be the ground state with a 
large basin of stability. 

The ground state for 0.8 < Ad < 1 appears to break 
the glide plane symmetry while remaining invariant un- 
der C2 rotations. While larger systems would be required 
to make definitive conclusions, this evidence is in stark 
contrast to the prediction that the kagome lattice ground 
state is a spin liquid^^'^". 

Remarkably, on this lattice all valence bond crys- 
tals break the glide plane symmetry (the most glide- 
plane symmetric configuration of Fig. [UJb) still breaks 
glide plane symmetry if the missing dimers are added to 
the picture). So the breaking of glide plane symmetry 
strongly supports recent proposal a^i^i^^ that the ground 



state maybe a valence bond crystal (VBC). Candidate 
such states include the Marston and Zeng 36 site unit 
cell dimer (spin singlet) pattern(MZ)i^ and the pin-wheel 
pattern (in the presence of distortion) (see Fig. [5]). One 
may argue that the MZ pattern should also break C2 
rotational symmetry. However, such a symmetry may 
naturally be restored by benzene-like resonances on the 
three dimer hexagons. One should note that whether the 
recent ED resultsi^ are disfavoring the MZ VBC state as 
well as other proposed VBC statesii^ or not, is a subject 
of intense debatei^^ 

Since any VBC will have a diverging glide-plane sus- 
ceptibility the results presented here are not very sen- 
sitive to transitions between different VBC's as long as 
they conserve C2 symmetry. For example, a transition 
from the MZ pattern at Ad = 1 to the pin- wheel pattern 
at Ad < 1 is certainly possible. For 0.8 < Ad we can 
therefore not exclude the presence of several different C2 
symmetric VBC phases although the rate of divergence 
of xgp could potentially be quite different for different 
phases. In fact, one might speculate that the cusp in xgp 
for = 36 and in XC2 for A^ — 24, 36 in both cases at 
A^ = 0.98 is a signature of a phase-transition between 
different valence bond crystals. 

A phase transition near Ad ~ 0.8 was also found in the 
large-N study of Ref. Thus both large-N and exact di- 
agonalization methods predict the existence of a quantum 
phases transition at a value of Ad away from the ideal 
kagome limit. If we assume that the spin gap is non-zero 
at Ad = 1 and vanish approximately linearly with the 
deviation of Ad from 1, then this value for the quantum 
critical point is also roughly consistent with the vanish- 
ing of the spin gap (which, from exact diagonalizations, 
is estimated to rather small but finit o^°'^^ in the ther- 
modynamic limit and has a value of 0.1848J at Ad = 1 
for the 36-site cluster). Given the apparent glide plane 
symmetry breaking for Ad ^ 0.8, this phase transition 
then appears to be between two phases with unrelated 
orders. It may therefore fall into the deconfined univer- 
sality classic if it is not a first order transition. 
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